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Abstract. - We study the effect of advection and noise on the field theory for directed percolation 
(DP). We show that even a very small advective velocity is enough to change the universality 
class of the dynamic phase transition. When the noise is taken to be proportional to the square 
root of the population density, we find an additional nonequilibrium "spinodal" line separating a 
region where an exponentially decreasing density is metastable, from another one in which it is 
unstable. If the noise is instead linear in the density, the phase diagram changes dramatically both 
quantitatively and qualitatively, and the spinodal line becomes a true phase boundary. We briefly 
discuss possible applications of our results to microbial sedimentation and population dynamics 
in rivers. 
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The concept of universality classes has been employed 
with much success in equilibrium statistical physics, gath- 
ering myriad phase transitions into a handful of classes. 
jThe effort to extend and apply this concept to nonequi- 
librium systems is still ongoing, and is one of the fore- 
most challenges of statistical physics. Probably the widest 
and best characterised nonequilibrium universality class to 
date is that of directed percolation (DP) [TJ, a system of- 
ten described by the following Langevin equation [2]: 



d t p = Dd xx p + ap- bp 2 + ^/prj 



(1) 



X 



Here, p — p(x,t) is a density, D a diffusion coefficient, 
a > a "growth rate", b > a saturation constant, and 
r)(x, t) is Gaussian white noise of unit variance. 

DP is an archetypal nonequilibrium phase transition, 
which takes an active phase, where the fluctuating den- 
sity is non-zero, continuously into an absorbing one with 
zero density. It succesfully describes systems as diverse 
as chemical reaction-diffusion processes, epidemic spread- 
ing, percolation through porous media, growing microbial 
population and branching-annihilating random walks [3]. 
This broadness is now expressed by the "DP conjecture" , 
stating that all systems exhibiting a continuous transi- 
tion into a unique absorbing state, characterised by a one- 
component order parameter, and not showing any extra 
symmetries or conservation laws, belong to the DP uni- 
versality class 0H]. The DP conjecture has proven ex- 
tremely robust to changes in the microscopic dynamics. 



This brings us to the subject of this article: what is the 
fate of the DP universality class in the presence of advec- 
tion? 

One might expect that the addition of an advection term 
to |T]) could be transformed away by a Galilean transfor- 
mation 0E] thus rendering it irrelevant. This is however 
not true in the presence of fixed boundaries. For most 
equilibrium systems boundaries are not expected to affect 
the critical behaviour, however in a nonequilibrium system 
boundaries are known to play a crucial role in determin- 
ing phase transitions [7J. To illustrate this and gain in- 
sight into our question, we consider the effect of advection 
on the noiseless limit of Eq. [TJ which is nothing but the 
celebrated Fisher-Kolmogoroff (F-KPP) equation 5 . It 
clearly admits two steady state solutions (p = 0, p = a/b) 
and exhibits the well-known Fisher Waves: wavefronts 
emerge from portions of the system in the high density 
state p = a/b and propagate into empty regions with ve- 
locity Vf = 2\jDa 5 . It has recently been shown 8,9 
that in the presence of boundaries, and upon addition 
of an advection term vd x p, F-KPP exhibits a new low- 
density steady-state: an exponential phase in which the 
density profile decays exponentially away from the bound- 
ary. The competition between the advancing Fisher wave 
and the advection term triggers a novel discontinuous non- 
equilibrium phase transition separating this low-density 
phase (v > Vf) from the high-density one (v < Vf), in 
which p ~ a/b throughout the system. Restoring noise 
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Figure 1: Phase boundary for a — 0.5 and g oc yj~p. The critical 
points found via the Dickman and Dornic et al. algorithms are 
shown as black crosses and red diamonds, respectively. Both 
methods agree up to an accuracy of about 1%. The solid ma- 
roon line correspond to the theoretical prediction (see below) 
for the spinodal line whereas the blue squares correspond to its 
numerical counterpart (see below). 

and returning to DP, it is natural to question the robust- 
ness of the exponential profile. Indeed, we now have two 
candidates for the low-density phase: the exponential one 
or the usual absorbing state of DP. 

In this work, we show that in the presence of advection, 
the low-density phase is always absorbed at long times 
whereas the fate of the exponential state depends on the 
strength of the noise; it is never seen at large noise but 
is metastable below a certain threshold. The role played 
by the noise in this dynamical transition is crucial. We 
explore this issue further by considering noise proportional 
to the density rather than to its square root and show that 
this transforms the dynamical transition into a true phase 
transition. Most importantly, we find that the addition of 
an advective term to DP, however small, always changes 
the nature of the transition, rendering it discontinuous. 
In other words, advection is a relevant perturbation in 
renormalization group jargon, in a similar way that, for 
instance, a magnetic field is a relevant perturbation in an 
Ising model in equilibrium statistical mechanics. 

Our starting point is the following extension of Eq. (JXJ) 

d t p = Dd xx p + vd x p + ap- bp 2 + T a g(p)ri (2) 

Our study is entirely carried out in one dimension; ex- 
tension to higher dimensions would be interesting but, 
since the advection term can be extended in various ways, 
we leave this for further study. We consider two differ- 
ent kinds of nonequilibrium noise distinguished by the 
functional form of g(p). In the first part of the article, 
we choose g(p) = -J p(x, t), as in conventional DP (see 
Eq. ([T])), and refer to this as "square root" noise. Such 
noise typically arises from fluctuations among individuals 
in a finite population of average density p{x, t). In the sec- 
ond case, that we address toward the end of the paper, we 
consider a "linear" noise, defined by g(p) = p{x,t), which 
is usually met in the context of population dynamics to 
model fluctuating environments [TU] ■ This kind of noise is 



sometimes called "multiplicative" in the literature [3] , but 
we use this terminology here to describe generic depen- 
dence of the noise on the density. Once g(p) is chosen, Tq 
is a parameter used to scale the strength of the noise rela- 
tive to the other terms. Throughout the paper, we impose 
'no-flux' boundary conditions at the top and bottom of 
the system^, whose size is L. As is usual for semi-infinite 
systems, we first take the large time limit before sending 
the system size to infinity. 

To obtain our numerical results for the square root 
noise, we used two different algorithms. The main techni- 
cal difficulty is to ensure that, upon time discretization of 
Eq. ©, the multiplicative noise does not lead to unphys- 
ical negative densities when p is small (which is the case, 
for instance, when we are close to the critical point of DP 
at v = 0). The first method we use was originally pro- 
posed by Dickman [11) , and entails a symmetrical trunca- 
tion of the noise coupled with discretization of the density 
such that TQ^fprj never exceeds —p(x,t). The second is 
an improved integration scheme proposed by Dornic et al. 
in |12j . Instead of drawing a random Gaussian number at 
each x to be later on multiplied by Toy/ p(x, t), this al- 
gorithm draws directly from the probability distribution 
function which solves the Fokker-Planck equation associ- 
ated with the linear part of the local Langevin equation at 
each x (see Ref. [H] for details). After this, the determin- 
istic non-linear part of the equation is evolved via a finite 
difference scheme. In this case, the density is never nega- 
tive by construction and numerical problems are avoided. 
We checked that both schemes yield very similar results. 
In all the simulations presented in this article, we use 
a = 0.5, 6=1 and D = 1. Typical integration param- 
eters with Dickman's algorithm are dx = 0.1, dt = 0.001, 
L = 1000 whilst varying v and Tq. For the Dornic et al. 
algorithm, whose generalisation to an advection term is 
quite straightforward, we could choose larger dx up to 0.5 
and dt up to 0.05, considerably reducing simulation time. 

First, we map out the phase diagram in the To — w pa- 
rameter plane (see Fig. [I]). A full line of critical points 
v = v c (Tq) now links the two limiting cases of DP (v = 0) 
and F-KPP with advection (Fo = 0). For small To and 
v, we observe a high-density phase, where the total mass 
M of the system is extensive with system size, whereas 
the rest of the phase diagram is composed of low-density 
regions (whether the exponential profile or the absorbing 
state) where M/L — > as L — > 00. The absorbing state 
can be accessed via fluctuations and is thus favoured by 
larger values of Tq. As a consequence, the critical ve- 
locity v c (Tq), above which the high density phase is not 
observed, is a decreasing function of noise strength. 

Interestingly, long time simulations of the system show 
that its dynamical behaviour in the low-density phase is 
not uniform. The exponential phase is long-lived for small 
Tq and we refer to this regime as the 'low noise' case. 
While holding the velocity fixed, increasing Tq has the 

i. e, we take D8 x p + vp = at x = 0, L 
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Figure 2: Plots of the average order parameter (m) as a 
function of v, for To = 0.10 and for three different system 
sizes (see legend). The presence of a crossing of the order 
parameter curves indicates that the transition is discontinu- 
ous. Inset: probability distribution of the order parameter, for 
v = v c ~ 1.21 and L = 100. 



Figure 3: Plots of the average order parameter (m) as a 
function of v, for To = 0.50 and for three different system 
sizes (see legend). The presence of a crossing of the order 
parameter curves indicates that the transition is discontinu- 
ous. Inset: probability distribution of the order parameter, for 
v = v c ~ 0.58 and L = 100. 



effect of decreasing the average lifetime of the exponential 
state, until we can barely see it; this is the beginning of 
the 'strong noise' case. As we shall see, this threshold of 
To corresponds to a (nonequilibrium) spinodal line which 
separates the two dynamical regimes (See Fig.fTJ). We now 
explain how the critical points were obtained and discuss 
the nature of the phase transition, detailing in particular 
the differences between the low and strong noise regimes. 

In order to identify the critical points we define the or- 
der parameter m — M/L, where M is the total mass in 
the system. We compute the average of m in the quasi- 
stationary stated over a number of different simulations, 
and then plot this average as a function of v, for fixed Tq. 
For finite systems, the crossover between low- and high- 
density phases sharpen when the system size is increased, 
but the order parameter curves intersect at a well-defined 
non-zero value of v (see Fig. 2 and 3), which we take as the 
critical velocity v c . This scenario is indicative of a discon- 
tinuous transition: in the thermodynamic limit, as v /*• v c 
the density in the stead state approaches a non-zero value, 
whereas for v \< v c it is strictly zero. Thus, for all non- 
zero values of v, there will be a discontinuous "jump" at 
the transition point; the system no longer belongs to the 
DP universality class. 

Whereas the nature of the transition differs from that 
of DP, it is also different from the one of the noiseless case 
(ro = 0) discussed in the introduction. In the absence 
of noise, a band structure is formed at criticality and a 
stationary front separates high- and low-density regions 
in the steady state. In the noisy case considered here, 
the transition is in general different. Even the low noise 
regime, while still bearing some resemblance to the noise- 
less case, shows qualitatively distinct behaviour. This last 
point can be appreciated by looking at the quasistation- 

2 This is the steady state of the probability distribution condi- 
tioned on survival. 



ary probability distribution of the average density of the 
system at criticality, which is shown in the inset of Fig. 
2. This distribution is composed of an approximately flat 
part coexisting with a peak close to zero. The physical 
interpretation is that in the low noise regime, the band 
is still present but the front, no longer stationary, now 
performs a random walk and occasionally gets stuck in 
the exponential profile. We checked numerically that the 
roaming band indeed performs a random walk by mea- 
suring the variance of the order parameter as a function 
of time (data not shown). The order parameter distribu- 
tion is linked, in equilibrium statistical mechanics, to the 
shape of the free energy. In the case of a standard discon- 
tinuous transition, close to criticality, there would be two 
coexisting free energy minima which result in two peaks 
in the probability distribution of the order parameter. In 
our case however, we have coexistence between one peak 
corresponding to the exponential profile and a flat piece in 
which all values of the total mass are essentially equiprob- 
able. The latter piece of the distribution corresponds to 
the high-density phase, in which the average of the order 
parameter thus scales with the system size. 

The transition into the low-density, large noise regime 
is quite different and, as might be expected, is reminiscent 
of what happens in DP. In this regime, there is no longer 
any front, or band structure, due to the large fluctuations, 
and the system now uniformly collapses into the absorb- 
ing state. The numerics in this region are more difficult. 
Since we are working with finite size systems, when v < v c 
the system is occasionally absorbed, even though it prop- 
erly belongs to the high-density phase. To overcome this 
difficulty, we use PERM (Pruning and Enrichment Rosen- 
bluth Method, see [T3] for a recent review) for simulating 
the quasistationary state. Despite the apparent similarity 
with DP, the transition is still discontinuous, as can be 
seen in Fig. 3. However, the order parameter distribu- 
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tion at criticality also differs from the low noise case (see 
the inset of Fig. 3). First, the peak corresponding to the 
exponential phase has disappeared. Then, the flat part 
in the distribution is now replaced by a broad peak at a 
finite value. For completeness, the quasistationary distri- 
bution, obtained using PERM, should be complemented 
by a delta function at (p) = 0. Therefore this case is closer 
to a standard discontinuous transition, with two compet- 
ing peaks at criticality. In the low density phase, only the 
delta function survives in the long-time limit. 

As previously mentioned, the frontier between the low 
and strong noise regimes can be probed via the dy- 
namical stability of the exponential phase. To do so, 
we study the evolution of the total mass in the system 
M(t) = Jq p(x,t)dx, obtained by integrating eq. @ with 

g(p) = Vp : 



M = aM-b / p 2 (x,t)dx+T Q / y/ p(x,t)r](x, t)dx (3) 
Jo Jo 

The diffusion and advection terms have dropped out due 
to the boundary conditions. We approximate this by the 
following simplified Langevin equation^: 



(4) 



M = aM - j3M 2 + r VM))(«) 



where f](t) is Gaussian white noise of unit variance. Both 
noise terms in ^ and (U]) are gaussian, have the same 
mean and variance, and are therefore equivalent. The 
term —j3M 2 is an approximation — we know from ((3]) that 
there must be saturation terms in the effective dynamics 
of M(t), and retain only the lowest order in M. We believe 
this approximation to be reasonable as we only consider 
the low-density regime, and hence small mass. (In the 
high-density phase, the mass would be extensive with the 
system size and the approximation would break down.) 
The parameter /3 contains the dependence on v and b but 
its form is not known exactly. 

The dynamical stability of the exponential is not easily 
studied from (j?]) as the noise is multiplicative. We there- 
fore recast it into an additive Langevin equation via the 
following change of variable u = 1\f~M /Y. Using the Ito 
formula [33], ^ becomes 



. . au 2 PY 2 u 4 

V cS (u) = h 

V ' 4 32 



loe 



u = -d u V eS (u) + rj; 

(5) 

The problem is now reduced to the diffusion of a particle 
in a potential V c s at temperature T = 1/2. A steady- 
state solution is thus given by P(u) oc exp[— 2V c ff(u)]. One 
notes however that P(u) ~ w -1 when u — > 0. The poten- 
tial is thus not normalizable, which simply stresses that 
the low density phases are always absorbed as t —¥ 00; the 
only normalizable steady-state solution is P(M) = S(M). 
The shape of V e s nevertheless contains relevant informa- 
tion for the dynamics, as illustrated in Fig. 2k- When 



ro < r c = -^p, the effective potential has a local min- 
imum corresponding to a potential well for positive M. 
The well lies above the global minimum at M — and 
corresponds to a metastable phase with finite mass: the 
exponential phase. Conversely, for Yq > Y c there is no 
metastable state and the system falls directly into the ab- 
sorbing state. T c thus corresponds to a spinodal point 
at which the exponential phase turns from metastable to 
unstable. 




Figure 4: (a) Plot of the effective potential for 77 oc ^/p, a = 2, 
/3 = 0.5. From top to bottom, we used Y = 5.0; 2.8; 2.0; 1.5 
while T c = 2V2. (b) Same for r\ oc p, a = 2, /3 = 0.5. From top 
to bottom, we used T = 0.5; 1.0; 2.0; 4.0 while V c = 2. 
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3 A similar approximation was used by Munoz in 3 



Figure 5: Contour lines in the (v, To) plane for which r equals, 
from left to right, 500, 100, 50 and 30. In order to plot the 
contour lines of the 0D model, we use j3 = 2.3 v (see below). 
The analytical prediction for the spinodal line is also plotted, 
and closely matches the curves for which r = 500. 

To construct the spinodal line r c (u) from the 0D model, 
we try to relate j3 to v (the dependence on b is not relevant 
since b is constant throughout our study). To do so, we 
compare the mean time r taken, for both the system and 
its 0D approximation, to reach the absorbing state. We 
compute numericalljO a set of contour lines for r in the 
planes (w,r ) and ((3,Y ), corresponding to the ID and 
0D systems respectively (see Fig. 5). Strikingly, simply 
setting (3 ~ 2.3 v suffices to make the two sets of con- 
tour lines overlap. This strongly supports the validity of 
the 0D model and hence validates our interpretation of 
the frontier between low and strong noise regimes as a dy- 
namical phase transition. We see on Fig. 5 that the the- 
oretical prediction for the spinodal line corresponds to a 

4 For the 0D model, this mean first passage problem can be solved 
exactly up to a numerical integration |15| 
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mean first passage time (MFPT) to absorption of t — 500. 
In simulations, we indeed observe long-lived exponential 
profiles when the MFPT to absorption is larger than 500 
whereas they are barely seen otherwise. The simulation 
data corresponding to r = 500 was thus used to pinpoint 
the spinodal on Fig. 1. 

We now turn to the investigation of the dependence of 
the phase diagram on the type of noise used in eq. ([2]), 
by considering g(p) — p. Without advection (v — 0), this 
model was studied in [3l|T6l [TT] and it was shown that 
the system can be mapped onto the KPZ equation, by 
means of a Cole-Hopf transformation. As a consequence, 
the absorbing state phase transition at v — is contin- 
uous, but does not belong to the DP universality class. 
For instance, holding Tq constant and varying a close to 
criticality yields (p) ~ (a — a c )^ with f3 = 1.5, whereas 
(3 ~ 0.22 in DP [TT]. In our system, upon switching from 
the square root to the linear noise, the phase diagram is 
dramatically altered (Fig. [B]). 

First, the exponential profile is now completely stable 
for small Tq and the transition between low and strong 
noise is a true phase transition. To understand why, we 
rely again on a 0D model, obtained by replacing \f~M by 
M in the noise term of eq. ((H). We this time consider 
the change of variable w = r _1 logAf to obtain an ad- 
ditive Langevin equation. Using the Ito formula [14] the 
corresponding equation reads 



w = -d w V eS (w) + r]; V eS {w) 



r2 < 



(6) 



and the putative steady-state is given by P(w) oc 
exp[— 2V e g (w)]. For T < \/2a, this is normalizable and the 
system is thus not absorbed, there is a proper normaliz- 
able steady-state distribution with M ^ 0. For T > \[2a, 
exp[— 2I4ff(iu)] is not normalizable and the steady-state 
distribution once again corresponds to a delta function, 
the system will be absorbed. Considering now the stability 
of the exponential profile, we see that the effective poten- 
tial has a single minimum which switches from w — — oo 
(M = 0) to finite w depending on the sign of the first term. 
The transition point r c = \/2a corresponds to the normal- 
ization criterion, the spinodal line has become a true phase 
transition. Note that the critical line r c (a,/3) = \/2a is 
now independent of /3 and hence of the velocity. Numeri- 
cally, we find that, for a = 1/2, the transition line is almost 
vertical and very close to T c (v) = 1. This close agree- 
ment with our predictions is reinforced by the absence of 
any fitting parameter. As we approach the high-density 
phase, the mass of the exponential profile increases, and 
the slight bend of the transition line is presumably due 
to non-linear terms beyond the quadratic one, neglected 
in the 0D approximation. In the low noise regime, for 
To < r c , the exponential state is now stable rather than 
metastable and hence constitutes a true phase. 

Beyond the transformation of the dynamical transition 
into a true phase transiton, the transition between the 
high- and low-density phases (whether exponential or ab- 




Figure 6: Phase diagram for r\ oc p (a = 0.5). The red squares 
correspond to the transition between high- and low-density 
phases. Within the low-density region, there is now a second 
true phase transition between low and strong noise, indicated 
by the maroon triangles. For comparison, we include the phase 
boundary of the square root noise case (black crosses). The 
horizontal magenta line corresponds to the zero-noise transi- 
tion point v c ~ 1.41 whereas the vertical one is the theoretical 
prediction from the 0D model. 



sorbing) is also changed significantly. To pinpoint this 
critical line, we proceed numerically as before. The model 
is now easier to simulate since close to the transition, when 
p <C a/b, the fluctuations remain of order p and are much 
smaller than in the previous case, where they scaled as 
yfp. Rare events leading to absorption are unlikely and 
Dickman's algorithm works very well, in particular dis- 
cretization of the density is unnecessary. Our results sug- 
gest that for small noise, fluctuations are irrelevant and 
one recovers a transition identical to one of the determin- 
istic limit — see Fig. [7] Indeed, numerically the transition 
line between exponential and high-density regimes is inde- 
pendent of To — hence horizontal — with v c ~ 1.41 as in the 
noiseless limit. Note that this transition is discontinuous 
and the band is stationary - in sharp contrast with the 
yfp noise where it performed a random walk at criticality. 
In the large noise case, where the stable state is absorb- 
ing rather than exponential, our simulations are consistent 
with a continuous transition (see Fig. 8) and more data 
would be needed to discriminate it from the v = limiting 
case. 

In conclusion, we have studied the effect of advection 
and the role of noise in the context of continuous absorb- 
ing phase transitions, and directed percolation in partic- 
ular. Our main result is that as soon as an advective 
velocity is added, the transition between the low- and 
high-density phases becomes discontinuous and therefore 
is no longer in the DP universality class. It is important 
to note that, in line with other terms in DP, the advec- 
tion term is a nonequilibrium term, as it is related to the 
transport rather than to the thermodynamic properties 
of the system. Furthermore, we have illustrated that the 
type of noise used to model non-equilibrium phase tran- 
sitions can control the kind of transitions observed. For 
instance, switching from square root to linear noise in our 
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Figure 7: Plots of the steady state density profiles for the linear 
noise (Fq = 0.1, solid lines) and noiseless (Fo = 0, dotted lines) 
cases, for L = 40 and values of v, from right to left, 1.30, 1.36, 
1.39, 1.40 and 1.41. 
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Figure 8: Plots of the average order parameter (m) as a func- 
tion of v, for To = 1.25 in the linear noise case, for three dif- 
ferent system sizes (see legend). This figure is consistent with 
a continuous transition. 



system transformed a dynamical phase transition into a 
true phase transition, stabilizing a previously metastable 
state. In addition to this, it appears that a large por- 
tion of the critical line has changed from discontinuous to 
continuous, although more extensive simulations would be 
required to confirm this. This important quantitative and 
qualitative difference highlights once again [3] that care 
should be exercised when deriving fluctuating hydrody- 
namic equations for non-equilibrium models, as the very 
form of the noise, which is sometimes overlooked, may 
drive unexpected changes in the physics of the system. 

We close with a brief mention of a possible application 
of our results to experiments with microbial systems. In 
population dynamics, Eq. [2] at zero or small noise has 
been used to model bacterial sedimentation in pipettes [5] , 
plankton sinking in oceans [18 (although in that work the 
logistic growth law was substituted by a non-local forcing 
term) , and the persistence of populations of aquatic organ- 
isms in rivers [19] (i.e. the problem of how organisms can 
resist being swept upstream by advection from the flows 
in the rivers they inhabit). In all of these cases, experi- 
ments suggest that the Fisher equation is a good, although 
approximate, starting point to study these systems in the 
absence of advection: therefore the regime of interest in 
our calculations is the small noise one where a front can 
still be identified. Our results suggest that in this situa- 
tion, being at or (for finite size systems) close to criticality 
brings with it a giant increase in fluctuations, and as a re- 
sult we would predict that the system under these condi- 
tions would be unable to reach a steady state: the Fisher 
wavefront would indefinitely perform a random walk in- 
stead. It would be interesting to see whether controlled 
experiments in bacterial colonies subject to, e.g., uniform 
flow in microfluidic devices may be designed which give a 
transition between a persistent phase in which the Fisher 
wave moves faster than the flow to another one in which 
the bacteria are swept away by the moving fluid. If this 
is the case, our results suggest that it would be extremely 
interesting to closely monitor the behaviour of these ex- 
periments close to this transition. 
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